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FIRST QUARTER SECOND YEAR PROGRESS REPORT 

This report covers technical progress during the first quarter of the second year of NASA 
Sun-Earth Connections Theory Program (SECTP) contract “The Structure and Dynamics of 
the Solar Corona and Inner Heliosphere,” NAS5-99188, between NASA and Science 
Applications International Corporation, and covers the period August 16, 2000 to December 15 
2000 Under this contract SAIC and the University of California, Irvine (UCI) have conducted 
research into theoretical modeling of active regions, the solar corona, and the inner heliosphere, 

using the MHD model. 

The renewal of Option 2 of this contract ( which renewed Year 2 of the contract from 
8/16/2000 until 8/15/2001 ), was delayed due to contracting delays at NASA. The renewal was 
not received until November 27, 2000. An a consequence, the progress reported in this 
« quarterly ” report includes a four-month period, rather than the customary thre e-month 
period. The actual funded time period during this time was from November 27, 2000 to 

December 15, 2000. 

In the following sections we summarize our progress during this reporting period, u 
descriptions of our work can be found in the cited publications, a few of which are attached to 

this report. 

Publication on the Magnetic Field Topology in Prominences 

The paper “Magnetic Field Topology and Modeling of Loops in Prominences,” by R. 
Lionello, Z. Mikic, J. A. Linker, and T. Amari, has been submitted for publication in The 
Astrophysical Journal Utters. We present a computational model of the magnetic field lines m 
a prominence. We describe how the process of flux dispersal and cancellation in a realistic 
configuration can lead to the formation of a stable flux rope structure. The flux rope fulfills 
several theoretical and observational constraints associated with prominences: twist, shear along 
the neutral line, and dips. We have performed ID hydrodynamic simulations along selected 
field lines in the final state to show that the flux rope can support cold and dense material t a is 
characteristic of a prominence. The referee’s comments have been received, and the paper is 
currently being revised for resubmission. This paper is included in Appendix A. 

Publication on an Empirically-Driven Global MHD Model of the Solar Corona and 
Inner Heliosphere 

The paper “An Empirically-Driven Global MHD Model of the Solar Corona and Inner 
Heliosphere,” by P. Riley, J. A. Linker, and Z. Mikic, has been submitted for publican on m 
Journal of Geophysical Research. We have developed a three-dimensional MHD model of the 
solar corona and heliosphere. We split the modeling region into two distinct parts, the so ar 
corona (1 solar radius, R s , to 30 R s ) and the inner heliosphere (30 R s to 5 AU). This combine 
model is driven solely by the observed line-of-sight photospheric magnetic field and can thus 
provide a realistic global picture of the corona and heliosphere for specific time periods of 
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interest. We have used the model to illustrate heliospheric structure during three different 
phases of the solar cycle: (1) Carrington Rotation (CR) 1913 (August 22, 1996-September 18, 
1996), which occurred near solar minimum and overlapped the “Whole Sun Mont 
campaign; (2) CR 1892 (January 27, 1995-February 23, 1995), which occurred during the 
declining phase of cycle 22 and coincided with the so-called “Ulysses rapid latitude scan”; and 
(3) CR 1947 (March 7, 1999-April 4, 1999), which occurred approximately 15 months before 
the predicted maximum of solar cycle 23. We compared Ulysses and WIND observations with 
the simulation for CR 19 13 and compared Ulysses observations during Us traversal from pole to 
pole with CR1892. We found that the simulations reproduce the overall large-scale features o 
the observations. Using the near solar maximum results, we have speculated on the structure of 
the high-latitude solar wind that Ulysses will encounter during its traversal of the southern and 
northern solar poles in 2000 and 2001 , respectively. In particular, the results suggested that, 
due to the presence of equatorial coronal holes, the ordered pattern of CIR tilts and their 
associated shocks, which was observed during Ulysses initial southward excursion in 1992, wi 
likely disappear completely as Ulysses moves toward the south pole. We anticipate that 
Ulysses will encounter fast streams but will not remain within them for more than a fraction of a 
solar rotation. Finally, the simulations suggested that crossings of the HCS will persist up to at 
least -70° heliographic latitude. The paper has been revised, based on the referee’s comments, 
and has been accepted for publication. This paper is included in Appendix B. 

Publication on Self-Consistent Formation of Prominences 

The paper “Magnetohydrodynamic Modeling of Prominence Formation Within a Helmet 
Streamer,” by J. A. Linker, R. Lionello, Z. Mikic, and T. Amari, has been submitted for 
publication in Journal of Geophysical Research. This paper describes the self-consistent 
simulation of the formation of prominences in the solar corona. We present a 2D axisymmetnc 
MHD model to self-consistently describe the formation of a stable prominence that supports 
cool, dense material in the lower corona. The upper chromosphere and transition region are 
included in the calculation. The prominence is formed by shearing a coronal arcade and 
reducing the magnetic flux, to form a magnetic field configuration with a flux-rope topology. 
The prominence forms when dense chromospheric material is brought up and condenses in the 
corona. The prominence sits at the base of a helmet streamer structure. The dense material is 
supported against gravity in the dips of the magnetic field lines in the flux rope Furt er 
reduction in magnetic flux leads to an eruption of the prominence, ejecting material into the 
solar wind. The paper has been revised, based on the referee’s comments, and has een 
accepted for publication. A copy of this paper is included in Appendix C. 

Eruption of a Magnetic Field Arcade by Flux Eruption/Cancellation 

It is well documented that strongly sheared magnetic fields observed on the photosphere 
are often associated with solar activity, some of which results in violent energy release in the 
form of coronal mass ejections and solar flares. The modeling of these events involves: (a) the 
buildup of the sheared field, (b) the triggering of the release of free energy, and (c) t e 
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subsequent dynamical events, including energetic plasma motions. We have demonstrated the 
energy buildup and release mechanisms by following the dynamic evolution of a magnetic 
arcade using 3-D numerical simulation. 

The sheared magnetic-field structures in the solar atmosphere can be dynamically formed 
in two fundamentally different ways as shown in our previous investigations. A current- 
carrying magnetic loop, for example, can be formed by plasma flows on the photosphere with 
an initially potential field (Van Hoven et al. 1995). The fluid motions induce an electric current 
by twisting the in situ field lines. Alternatively, a current-carrying loop can emerge directly 
from underneath the photosphere by breaking through the surface (Mok et al. 1997). The end 
results are remarkably similar, provided the physical parameters are comparable. In the present 
effort, we chose to build up the sheared magnetic structure by prescribing plasma flows on the 
surface at the location where there is initially a potential magnetic arcade field. The flow pattern, 
in general, can be arbitrary on the solar surface. However, we have identified two parameter 
regimes that lead to fundamentally different consequences in our dynamic model. In the first 
case, in which the spatial scale of the flow around the magnetic neutral line is small compared to 
the spatial scale of the field, free magnetic energy is built up, and stored, primarily in a small 
region along the neutral line. In the second case, the flows and the field have comparable spatial 
scales everywhere, and free magnetic energy is built up throughout the structure. The triggering 
mechanism to release the magnetic energy in our model is the emergence of a new, potential 
arcade field with polarity opposite to the original arcade. This process essentially weakens, and 
disrupts, the potential component of the sheared field structure. Since the plasma flows in the 
buildup phase are very slow compared to the Alfven speed, the sheared arcade is near force-free 
MHD equilibrium. The disruption of the potential component could trigger an analog of the 
ideal MHD kink instability , one of the most unstable modes often seen m laboratory devices. 

In the first parameter regime, we began with two elongated, parallel magnetic poles of 
opposite polarity as a potential field structure. A flux-preserving, sheared flow is applied to the 
surface in such a way that the flow is highly concentrated along the magnetic neutral line with a 
spatial scale much smaller than the dimension of the magnetic poles. The maximum flow speed 
at the base is 0.01 of the Alfven velocity, and this buildup phase lasts 190 Alfven transit times. 
At the end of this period, when the flows stop, free magnetic energy is built up primarily along 
the neutral line at a relatively low altitude in a highly localized manner. After the buildup, a new 
bipole with polarities opposite to the original is slowly emerged from the surface on a time scale 
of 100 Alfvdn times. The dynamic evolution is shown in Figure 1 . At r = 50 (measured after 
the beginning of the new flux emergence), the magnetic structure continues to nse slowly into 
the corona. Some of the field lines become helical as a result of reconnection between the new 
flux and the existing, strongly sheared flux, whose field lines are nearly parallel to the neutral 
line. The instability is triggered at approximately t = 80, when the strength of the new flux 
passes a threshold. The structure suddenly explodes upward on an ideal MHD time scale. The 
similarity of this eruptive behavior with the MHD kink instability is quite evident. By t = 87, 
the end of this simulation, the flux is being ejected out of the corona at a speed exceeding 0.4 
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times the Alfven velocity. This case is extremely dynamic. We believe that this is a plausible 
model for impulsive events such as coronal mass ejections and solar flares. 

In the second parameter regime, we applied a surface flow to the same bipole region to 
build up magnetic energy. The plasma flows, however, follow the contours of the vertical 
component of the magnetic field so that the field lines of the entire region are evenly subjected 
to the twisting motions. Free magnetic energy is built up throughout the structure, which 
expands slowly into the corona during the process. We were able to store more free energy in 
the system at the same value of resistivity, which eventually limits the growth of the electric 
current density. However, when the shearing stops and new flux is emerged, as in the last case, 
a magnetic tearing mode is excited on the vertical current sheet directly above the neutral line, 
resulting in an upward-downward flow pattern away from the X-point. The speed is typically a 
few percent of the Alfven velocity, and the evolution is on a slower time scale than in the 
previous case. This event is not explosive, and might correspond to the small plasma flows 
observed in some prominences. 

A preliminary version of this work was presented in the paper “Effects on Magnetic 
Structures by Disrupting Plasma Flows and Surface Magnetic Fields,’ by Y. Mok, J. A. Linker, 
and Z. Mikic, at the meeting of the American Astronomical Society/Solar Physics Division, held 
in Lake Tahoe, Nevada, June 18-22, 2000. These results are currently being prepared for 

publication. 

references 

Van Hoven, G., Mok, Y., & Mikic, Z. 1995, Ap. J ., 440, L105. 

Mok, Y., Van Hoven, G., & Mikic, Z. 1997 , Ap. J ., 490, L107. 
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Eruption of an Arcade Caused by Flux Emergence 



(a) / = 50 < b ) ' = 70 




(e) r = 85 



(d) i = 83 



Figure 1. A time sequence showing the dynamic evolution of the magnetic field lines of a 
highly sheared bipole region when new flux of opposite polanty emerges from underneath the 
photosphere. The new flux begins to emerge at t = 0. The sequence corresponds to (a) t 5 
thf r = 70 (cl t = 80 (d) t = 83, (e) t = 85 and (f) t = 87. The simulation was run on a mesh of 
127x79x71 points in (x,yz), with a Lundquist number of S = 5x 10 . The field lines in all six 
frames are traced from the same set of footpoints at the base. 
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APPENDIX A 


“Magnetic Field Topology and Modeling of Loops in Prominences 
R. Lionello, Z. Mikic, J. A. Linker, and T. Amari 
Submitted for publication in The Astrophysical Journal Letters 
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ABSTRACT 

We present a computational model of the magnetic field lines in a promi- 
nence, We show how the process of flux dispersal and cancellation in a realistic 
configuration can lead to the formation of a stable flux rope structure. The 
flux rope fulfills several theoretical and observational constraints associated with 
prominences: twist, shear along the neutral line, and dips. We perform ID hy- 
drodynamic simulations along selected field lines in the final state to show that 
the flux rope can support cold and dense material that is characteristic of a 
prominence. 

Subject headings: MHD — - Sun: corona - Sun: magnetic fields 
Sun:prominences 

1. Introduction 

Prominences are formed of material that is hundreds of times cooler and denser than 
the surrounding corona. They are often observed to lie along the neutral line separating 
regions of opposite polarity (Priest 1982). Considerable research has been performed in the 
past in order to find a favorable background magnetic configuration to support prominences 
against gravity (Priest 1982). One possibility is to start from a field generated by a dipole 
and apply a shear flow concentrated along the neutral line, which results in an equilibrium 
with dipped field lines (Antiochos et al. 1994). The presence of localized heating near the 
chromosphere makes it possible to form condensations, whose collective appearance is the 
prominence (Dahlburg et al. 1998; Antiochos et al. 1999). 

Our approach is driven by different theoretical considerations, which indicate that a 
flux rope is a viable candidate to support the denser and cooler prominence material. A 
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flux rope fulfills several of the theoretical and observational requirements such as twist, 
shear along the neutral line, and dips [however there are also observations that suggest 
the presence of vertical magnetic field in apparent contrast to this model (Zircher et aL 
1998)]. Amari et al. (1999) have shown how flux dispersal and cancellation in an idealized 
sheared arcade field can lead to the formation of a stable flux rope structure. In this work 
we show how a similar process leads to the formation of a stable flux rope and we study the 
condensation of plasma along the loop’s magnetic field lines in this topology. For our model 
we are motivated by observations of a prominence seen on September 23, 199G (see Fig. 1). 
The final configuration consists of a twisted magnetic flux tube embedded in an overlying, 
almost potential, arcade. We then solve the ID hydrodynamic equations along the field 
lines in this final state in order to show that the flux rope can support the cold and dense 
material that is characteristic of a prominence. Our work does not directly address how 
condensations are actually created. One possibility is that formation occurs by lifting of 
chromospheric material during the dispersal and cancellation phase. Siphon flows repiesent 
another possible mechanism (Antioc.hos et al. 1999) 

2. The MHD Model and Simulations 

Our approach is based on two separate models. We use a 3D MHD model to obtain 
the magnetic structure. Then we solve the ID hydrodynamic equations along field lines 
to obtain the plasma properties and show how the magnetic structure can support a 
prominence. Here we describe the 3D MHD model. 
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2.1. The MHD Model 

To model magnetic structures in active regions we solve the following MHD equations 
in Cartesian geometry: 


V X B = 

c 

(1) 

V x E = 

1 3B 
c dt 

(2) 

v x B 

E + = 

c 

r/J, 

(3) 

- 

'jxB + V' (f'pVv). 
c 

(4) 


E and B are respectively the electric and magnetic fields, v is the plasma velocity, c. the 
speed of light, p the resistivity, v the viscosity, p is the plasma density assumed to be 
constant. The exclusion of the pressure gradient and gravity from the momentum equation 
is justified because the dynamics of an active region is dominated by the strong magnetic 
field, i.e. the corona is a low-/? medium (where 3 is the ratio of the plasma and magnetic 
pressures). We use a 51 x 51 x 61 nonunifonn grid to mesh a cubic domain, -21 < x < 2 L, 
-2L < y < 2L, 0 < 2 < 4L, with L = 1.5 x 10 5 Km. The Lundquist number S is defined 
as the ratio of the resistive diffusion time r R to the Alfven time t a . In our model we have 
S = l x 10 3 , which is much lower than the value in the solar corona. This is necessary 
to dissipate structures that cannot be resolved since they are smaller than the cell size. 
However it does not influence the large scale dynamics in which we are interested. For the 
same reason we choose v such that the ratio of the viscous dissipative time \s. the Alfven 


time is = 200. 



2.2. The Simulation 


The MHD simulation can he divided into three main phases. In the first phase of our 
simulation we specify a magnetic flux distribution at the base of our domain. It consists 
of the superposition of three Gaussian distributions and captures the salient features of 
the bipolar region surrounding the prominence of September 23, 1996 (Fig. 1). We do not 
attempt to model this region in detail. We are merely seeking to model the large-scale 
magnetic topology as simply as possible. Then we find the corresponding potential magnetic 
field for this flux distribution, some representative field lines of which are shown in big. 2a. 

In the second phase we apply a velocity shear on the lower boundary surface for a 
specified duration, localized along the neutral line and advance Eqs. (3-4) in time. The 
form of the shear is 

v fc = V x v>(r, y)z, ( 5 ) 

where ib is a specified scalar function. If we choose ip — ip{B z ), it can be shown that the 
flow derived from Eq. (5) preserves the original magnetic flux of the configuration (Amari 
et at. 1996). The corresponding flow consists of two vortices rotating in the same direction 
and centered on the opposite polarity regions. The maximum speed of the flow is set to 
one hundredth of the Alfven speed. The system evolves through a series of force-free, 
increasingly sheared equilibria. With this flow pattern we do not intend to model a physical 
process by which the prominence is actually formed but we only use it as a method to form 
such equilibria with a given surface flux distribution. Figure 2b shows some field lines at 
the end of the shearing phase. 

Finally, in the third phase we reduce the magnetic flux at the photospheric boundary 
as in Amari et al. (1999). The magnitude of the original magnetic flux decreases by 20%. 
Our experience tells us that the exact change in the flux distribution is not important as 
long as the total flux is reduced. It is during this last phase that magnetic reconnection 
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occurs and a flux rope is formed (Fig. 2c). The flux rope, which is surrounded by an almost 
potential arcade, lies along the neutral line. The magnetic field in the flux rope has dips; wo 
will show how the cold and dense material tvpiral of a prominence can be stored in there. 


3. Dips and Condensations 


Our philosophy is to separate t he calculation in two parts. In the previous section wc 
have obtained the magnetic field with the 3D MHD model. In order to see if the magnetic 
field can support a condensation, wc now proceed to find the plasma properties, p and 7 . 
along selected field lines with a ID hydrodynamic model. This can be done because the 
plasma is low 3 and the perturbation to B due to pressure gradients and gravity is small. 


We solve the following equations: 


dp 

at 

ar 

In 


dp 

"ds 

ar 

ds 


dv dv 

ai + v as 


' a , 
-<’Aa; Av ' 


= -h - - <7 - 1)£ + --« T > - ") 


( 6 ) 

(7) 


1 dp 
p ds pA 


1 d ( dv\ 


( 8 ) 


s is the length element along the field line; T, p and p are the plasma temperature, pressure 
and density; k is Boltzmann constant; A(s) = Bq/B{s) is the expansion factor of the field 
line; ko is Spitzer coefficient for thermal conduction (ten = n 0 l j/2 ); 7 = 5/3 is the ratio of 
specific heats; Q is the radiation law [see Athay (1986)]; n„ and n p are the electron and 
proton number density and are equal in our calculation; H = 7/ 0 exp(-(r — /?^)/A) is a 
parameterized heating function; g, - g ■ b is the projected gravitational force; m p is the 
proton mass and p ~ rripTij,. 


According to Rosner et al. (1978), in the case of low-lying (constant pressure) loops 
there exists a lower limit L„ for the length scale A of an exponentially decaying heating 
function in order to have a solution with maximum temperature in the middle. This limit 
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is proportional to the length of the loop L. The presence of gravity does not fundamentally 
change this result (Scrio et al. 1981). When A < Ln only a solutioii with a temperature 
minimum in the middle is possible. The presence of a condensation in the middle can 
be thought as splitting the loop in two of length L/2. This effectively decreases Ln to 
Lu / 2 and if A > L u /'l then the solution with condensation in the middle is stable. Our 
choice for all simulations is A = 0.05 H ., . H 0 is chosen to give a surface energy flux of 
10 fi erg cur 2 s' 1 . 

As an initial state we choose an arbitrary temperature profile and calculate the 
corresponding pressure and density from solving Eq. (8) with v — 0. The temperature 
profile is chosen such that the system is not in thermal equilibrium. Equations (6-8) sire 
advanced to steady state using rigid wall boundary conditions at the ends of the loops (we 
have tested the stability of the final states by continuing the simulations with free-flow 
conditions (Antiochos et al. 1985; Mok ct al. 1991). With these boundary conditions at 
s = 0 and s = ia dense and cold chromosphere is formed, next to a transition region with 
sharp gradients. In the interior we have a hot, tenuous corona. Two kind of solutions are 
possible: with and without condensation. Loops with a height profile as in Fig. 3 (i.c., loops 
m the overlying arcade) do not support a condensation. In fact an eventual condensation 
forming at the top of the loop would be unstable to the Rayleigh-Taylor instability (dense 
material on top of tenuous). The solution is the classical loop solution with temperature 
maximum and density minimum at the top. These loops constitute the arcade surrounding 
the flux rope. However, when a dip is present in the magnetic field line as in Fig. 4, then 
the configuration is able to sustain the condensation. This occurs for the field lines that 
belong to the flux rope. 

A composite plot of all the simulations performed is presented in Fig. 5. The 
temperature and logarithmic density along the field lines is indicated in a color scale. 
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Condensations appear to be associated with the dips of the flux ropes and are absent in the 
semi-circular loops. The images can be compared with Fig. 1, ami they show how it is the 
superposition of several single condensations that forms the prominence. 

4. Conclusions 

We have shown how the process of canceling and dispersal of magnetic flux in a 
sheared realistic magnetic configuration can lead to the formation of a flux rope. Our 
simulation is intended to model the prominence that w'as visible on September 23, 1996. 
Solutions of the hydrodynamic equations with thermal conduction, heating and radiation 
losses were performed along the field lines of the final state. They show that the dips in the 
magnetic field lines in the flux rope can support the dense and cold material typical of a 
prominence. Arcade-like field lines do not have dips and are therefore incapable of confining 
the heavy material against gravity. The superposed view of all the condensations appears 
as a prominence. 

Real filaments have barbs also, which are curved feet connecting the body of the 
filament to the chromosphere. We plan to use our model with a more realistic magnetic 
flux distribution that includes parasitic polarities, which may be able to reproduce barbs. 
We have not addressed the point of how the prominence is actually formed, either via 
evaporation and successive condensation of chromospheric material or through mechanical 
lifting of cold and dense matter during the reconnection phase, ln order to investigate this 
issue a more sophisticated model is needed, which should include a self-consistent 3D MHD 
simulation with an energy equation incorporating coronal heating, thermal conduction and 
radiation losses. This will be studied in the future. 
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Fig. 1. Ha image of the prominence of September 23, 1996 from the National Solar Ob- 
servatory /Sacramento Peak, Sunspot NM (left). Corresponding magnetogram from the Kitt 
Peak National Observatory AURA/ National Solar Observatory /National Science Foundation 
(right). 

Fig. 2. - Magnetic field lines during the three phases of the simulation: (a) potential field; 
(b) localized shearing along the neutral line; (c) emergence of opposite polarity flux and flux 
rope formation, surrounded by an arcade. Underlying is the magnetic flux distribution. 

Fig. 3. — Height (top), temperature (middle) and density (bottom) profile as a function of 
the distance along the field line for a configuration that does not develop a condensation. 
No dip is present to sustain a stable condensation. 

Fig. 4. Height (top), temperature (middle) and density (bottom) profile a function of 
the distance along the field line for a configuration that develops a condensation. The 
condensation appears in the dip of the field line. 

Fig. 5. Temperature (top row) and density (bottom row) along several field lines and three 
different views from left to right. Arcade-like loops do not develop a condensation. Field 
lines forming the flux rope show that condensations occur in the dips. The combined view 
forms the prominence. 
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Abstract. 

In this study we describe a three-dimensional MHD model of the solar corona and 
heliosphere. We split the modeling region into two distinct parts: the solar corona (1 
solar radius. /?,, to 30 /?,) and the inner heliosphere (30 R, to 5 At ). This combined 
model is driven solely by the observed line-of-sight photospheric magnetic field and can 
thus provide a realistic global picture of the corona and heliosphere for specific time 
periods of interest. We use the model to illustrate heliospheric structure during three 
different phases of the solar cycle: (1) Carrington Rotation (CR) 1913 (August 22. 

1 996 - September 18, 1996). which occurred near solar minimum and overlapped the 
‘‘Whole Sun Month" campaign: (2) CR 1892 (January 27. 1995 - February 23, 199o). 
which occurred during the declining phase of cycle 22 and coincided with the so-called 
“Ulysses rapid latitude scan": and (3) CR 1947 (March 7, 1999 - April 4, 1999). which ‘ 
occurred approximately 15 months before the predicted maximum of solar cycle 23. We 
compare Ulysses and WIND observations with the simulation for CR 1913 and compare 
Ulysses observations during its traversal from pole to pole with CR 1892. We find 
that the simulations reproduce the overall large-scale features of the observations. We 
use the near solar maximum results to speculate on the structure of the high-latitude 
solar wind that Ulysses will encounter during its traversal of the southern and northern 
solar poles in 2000 and 2001. respectively. In particular, the results suggest that, due 
to the presence of equatorial coronal holes, the ordered pattern of CIR tilts and their 
associated shocks, that was observed during Ulysses initial southward excursion in 1992. 
will likely disappear completely as Ulysses moves toward the south pole. We anticipate 
that Ulysses will encounter fast streams but will not remain within them for more than 
a fraction of a solar rotaton. Finally, the simulations suggest that crossings of the HCS 
will persist up to at least ~ 70’ heliographic latitude. 
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1. Introduction 

Beyond ~ 10 solar radii (7t s ). solar material streams away from the Sun along 
roughly radial trajectories. A combination of temporal variations at the Sun. together 
with the rotation of the Sun. leads to parcels of plasma with different plasma and 
magnetic properties becoming radially aligned: faster material overtaking slower 
material leads to a compression front, while slower material being outrun by faster 
material leads to a rarefaction region, or expansion wave [Sarabhai. 19G3j. When the 
flow pattern at the Sun does not vary considerably during a solar rotation (such as 
during solar minimum), the large-scale compressive structures created by the interaction 
of these streams are fixed in a frame corotating with the Sun, and they are known as 
Corotating Interaction Regions (CIRs) [Smith and Wolfe, 1976], When the difference in 
speed between the slow and fast streams is sufficiently large, a pair of shocks may form 
bounding the CIR [e.g.. Pizzo. 1985). 

During its initial voyage out of the ecliptic plane. Ulyssses observed a systematic 
tilt to CIRs [Gosling el at. 1993. 1995a: Riley et al. , 1996]. These interaction regions 
were observed from January, 1992 to October. 1999. corresponding to the declining 
phase of solar cycle 22 In the a/imuthal (0) plane, their outward normals were tilted 
westward (i.e.. into the direction of planetary motion). This was to be expected, as the 
interaction regions aligned themselves roughly along the direction of the Parker spiral. 
In the meridional plane (0). their outward normals were tilted toward the heliographic 
equator. These results were supported by numerical MHD simulations [ Pizzo and 
Gosling. 1994] which suggested that these orientations were a natural consequence of 
the tilt of the solar magnetic-dipole axis relative to the solar rotation axis. This pattern 
of tilts was once again observed during Ulysses so-called ‘rapid latitude scan [Gosling 
et al.. 1995b] and most recently, as Ulysses descended from the nor then polar regions 
to the equator (during solar minimum and the earlv/mid ascending phase of cycle 23) 
[Gosling ct al., 1997: McComas et al., 1998a], 


Previous solar wind models have yielded considerable insight into the dynamical 
processes that shape the structure of the heliosphere. The largi-scale. time stationary 
structure of the heliosphere under an idealized tilted dipole geometry was studied by 
Pizzo [1991. 1994a]. At an inner boundary of 32 /?,. a flow pattern was constructed 
consisting of slow, dense flow about a heliomagrietic equator that was tilted by a 
specified amount relative to the rotational axis. Similar simulations were used by Gosling 
and Pizzo [1994] to interpret Ulysses observations at mid-heliographic latitudes, where 
it was found that CIR-associated reverse shocks persist to much higher latitudes than 
CTR- associated forward shocks. In both cases, the tilted dipole geometry at the inner 
boundary strongly drives the position and orientation of the CIRs and their associated 
shocks. In particular, under moderately low solar activity conditions. CIRs and shock 
outward normals are tilted toward the equator in both hemispheres. Thus, in the 
ecliptic plane, a spacecraft samples alternating tilts as it passes through a predominantly 
northern hemispheric CIR, and subsequently a predominantly southern- hemispheric 
CIR. Moreover, the strongest interactions take place away from the ecliptic plane, at 
latitudes roughly equal to the tilt of the solar dipole relative to the rotation axis. 

Wang and She.elcy [1990] exploited an empirical relationship between solar wind 
speed and coronal flux-tube expansion [ Levine et al.. 1977] to predict solar wind speed 
at 1 AU and beyond. They found that coronal flux tubes that expand more slowly 
correlate with faster asymptotic speed along that flux tube. Using this relationship, 
they predicted the types of wind speed patterns that Ulysses would be expected to see 
during its second solar orbit [Wang and Shecley, 1997] . In particular, they predicted 
that the ascending phase of the solar cycle would be dominated by recurrent (28-29 
day periodicity), high-speed streams originating from high-latitude extensions of polar 
coronal holes. Approaching solar maximum, however, persistent high-speed streams 
would disappear, only to be replaced by low -speed wind at all latitudes. Finally, at 
solar maximum (or more specifically at the time that the polar field reverses), very fast 
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episodic polar jets would be generated as active region fields migrated toward the solar 
poles. 

Usmanov et al. [2000] developed ail axisvinmetric MHD model of the corona 
and inner heliosphere. Their technique relies on outwardly propagating Alfven waves 
to provide all of the heating and acceleration of the solar wind flow. In reality, it 
is likely that a complex combination of thermodynamic processes is necessary to 
accurately describe the heating and acceleration of the solar wind [Lwnello et al., 2000]. 
Nevertheless, their approach quantitatively reproduces Ulysses' plasma and magnetic- 
field observations during the rapid latitude scan. Since the model is two-dimensional, 
no interaction regions are produced. In spite of this, their success suggests that such 
an approach in three dimensions might represent an effective compromise between a 
complex treatment that includes a transition region and coronal heating mechanisms 
[e.g.. Lwnello et al.. 2000] and a more limited polytropic approach [e.g.. Riley et al.. 
2000 ]. 

The heliospheric current sheet (HCS) is a global structure in the heliosphere that 
separates inwardly- and outwardly -directed magnetic field lines. It plays an important 
role in the modulation of cosmic rays and energetic particles in the solar wind. Simple 
models of the HCS have been constructed that ballistic-ally trace out the inferred 
locus of the neutral line at the Sun [e.g.. Joktpii and Thomas, 1981]. However, by not 
accounting for the deformation of the HCS through the evolution of interaction regions, 
these models can at best provide only a qualitative picture of the shape of the HCS. 

In contrast. Pizzo [1994b] used an idealized 3-D. steady-state MHD model to simulate 
the structure of the HCS out to 30 AU for the tilted-dipole geometries discussed above. 
He found that for dipole tilts larger than ~ 10°. noticeable deformation of the HCS 
occurred within 5 AU. and by ~ 10 AU (where CIR -associated shocks had overtaken 
the HCS) significant distortion of the HCS occurred. Pizzo emphasized the importance 
of the HCS as a global structure about which the interaction regions are organized. 


The data utilized in this study derives from the Solar Wind Over the Poles of the 
Sun (SWOOPS) ion sensor [Borne et al., 1992] and the magnetic field investigation 
[Balogh et al.. 1992] on board the Ulysses spacecraft and from the Faraday cup 
instrument on board WIND [ ( ) gih’ic et al.. 1995]. The plasma moments pioduced 
from the Ulvsses measurements have a typical resolution of 4-8 mins while moments 
produced from WIND measurements have a resolution of 90 secs. All data presented 
here, however, have been averaged over at least 1 hour. 

The purpose of this report is two-fold. First to introduce an empirically-driven, 
time-dependent. 3D resistive MHD model of the corona and heliosphere. To test its 
applicability, we make comparisons with WIND and Ulysses observations during a 
relatively quiet period of the solar cycle. We use the model to explore how the large-scale 
structure of the heliosphere changes over the course of a solar cycle, We then use these 
results to make basic predictions about the large-scale structure Ulysses may find at 
high heliographic latitudes during its second orbit at solar maximum. This report is 
structured as follows. In section 2 we introduce the basic features and limitations of 
the model. In sections 3. 4. and 5. we model the structure of the heliosphere at three 
time periods corresponding to solar minimum, the declining phase, and solar maximum 
conditions, respectively. In section 6 we explore the variability of the HCS for the 
aforementioned time periods. Finally, in section / we summarize the results of this 
study and speculate on what Ulysses may see when it traverses the poles of the Sun at 
solar maximum. 

2. Description of the Model 

2.1. Equations and Solution Technique 

In recent years, we have developed a 3-D, time-dependent resistive MHD model to 
investigate the structure of the solar corona [e.g.. Mikic et al., 1999: Linker et al., 1999]. 



8 


We solve the following system of partial differential equations, in spherical coordinates: 
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where B is the magnetic field intensity: J is the electric current density: E is the electric 
field: v is thp plasma velocity: p is the plasma inass density: p is the gas pressure: g is 
the acceleration due to gravity; 7 is the ratio of specific heats: i) is the plasma resistivity: 
u is the kinematic viscosity, and S represents energy source terms. 

In this study, we approximate the energy equation with a simple adiabatic energy 
equation (i.e.. S — 0) and choose the polytropic index (7) to be 1.05 in the coronal model 
and 1.5 in the heliospheric model The coronal value reflects the fact that temperature 
does not vary significantly in the corona. While this approximation significantly 
simplifies the calculation and reduces the time necessary to complete a simulation, the 
resultant plasma parameters predicted by the model do not show the same degree of 
variation as is inferred from in situ and solar observations. In the heliosphere, on the 
other hand, a polvtropic relationship between pressure and density is observed to hold, 
with 7 ~ 1.5 for protons (e.g.. Totten et at. [1995]: Feldman et al. [1998]). 

The details of the algorithm used to advance the MHD equations are provided 
elsewhere [Mikic and Linker. 1994: and Lionrllo et al.. 19981. Here we briefly make a 


few remarks. In the radial (r) and meridional (0) directions we use a finite-difference 
approach. In azimuth (<7). the derivatives are calculated pseudo-spectraily. We impose 
staggered meshes in r and 0 that have the effect of preserving V - B — 0 to within 
round-off errors for the duration of the simulation. The simulations discussed here were 
performed on {r.9.4>) grids ranging from 81 x 81 x G4 to 121 x 121 x 128. 

It is convenient to separate the region between the solar photosphere and the 
Earth into two parts. We distinguish between the more complex ‘coronal region, 
which includes the region from the photosphere up to 20-30 /?,. and the heliospheric 
region, which covers the region of space between 30 and 1-5 AL In the lattei. the 
flow is everywhere supersonic, gravity can be neglected (although it is included for 
completeness), and the energy equation can Ire reasonably approximated by a polvtropic 
relationship. Thus the time step required to advance the solution in the heliospheric 
model is considerably larger than the time step required for the coronal solution. 
Computationally then, it is more efficient to advance the heliospheric portion of the 
simulation independently of the coronal time step. 

Ultimately, our goal is to use the output of the coronal solution directly to provide 
the inner boundary condition of the heliospheric model. At present, liowev er. the 
plasma speed predicted by the polvtropic coronal model is not sufficiently high to 
drive the heliospheric model. The recent addition of a more complex thermodynamic 
treatment of the energy equation [ Lionello et al.. 2000] promises to yield more 
realistic plasma parameters: however, this refinement is still under development and 
requires further validation. Thus both as an interim solution, and as a practical and 
computationally-inexpensive approach, we have developed a technique for deducing 
speed, density, and temperature, based on the magnetic topology of the coronal solution 
(which is probably the most accurately determined property of the coronal model). 
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2.2. Modeling the Solar Corona 

Coronal calculations are typically performed between 1 /? 5 (i.e.. the base of the 
corona) and 20-30 /f 5 . although we have positioned the outer boundary as far away as 
1-2 AC [e.g.. Linker and Mikic. 1997). At the lower boundary, we specify the radial 
component of the magnetic held. B r . based on the observed line-of-sight measurements 
of the photospheric magnetic held, and uniform, characteristic values for the plasma 
density and temperature. An initial estimate of the field and plasma parameters are 
found from a potential field model and a Parker transonic solar wind solution [Parker, 
1903!. respectively. This initial solution is advanced in time until a steady state is 
achieved. 


2.3. Modeling the Heliosphere 

The heliospheric calculation is initiated using the results of the coronal solution. 
To determine the inner heliospheric boundary conditions (at r - 30/f s ). we use the 
magnetic field topology of the coronal solution. We assume that within coronal holes 
(i.e.. away from the boundary between open and closed magnetic field lines) the flow 
is fast. At the boundary between open and closed field lines, the flow is slow. Over a 
relatively short distance, we smoothly raise the flow speed to match the fast coronal 
hole flow. We specify the photospheric flow field as follows. 


iv(rf) = \ ~ + tiUlh ) J ( ' 

where d is the minimum distance from an open/dosed boundary, measured along the 
surface of the photosphere, a is a measure of how thick the slow flow band is (~ 0.1 


radians, or ~ O’), and w is the width over which the flow is raised to coronal hole values 


(~ 0.05 radians, or ~ 3°). Although this approach is somewhat ad hoc. it is based on 
the commonly held view that slow flow originates from the boundary between open 
and closed field lines and fast How originates in coronal holes. Once this speed map is 
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determined in the photosphere, it is mapped outward along field lines to generate the 
inner boundary of the heliospheric model at 30 R s . We emphasize that this approach 
is only a convenient method of generating the inner boundary conditions for the 
heliospheric solution. In reality, there is significant acceleration of the plasma between 
the photosphere and upper corona. 

In addition to the approximations discussed above, we make the following 
assumptions. First, we neglect the effect of pickup ions, which are thought to dominate 
the internal energy of the solar wind beyond C-10 AU [ Axford , 1972], Thus we limit our 
modeling region to < 5 AU. Second, we neglect the effects of differential rotation, which 
may play a role in connecting high latitude field lines near the Sun with lower-latitude 
interaction regions much further away ( Fisk, 1996], We assume that the inner boundary 
rotates rigidly with a period of 25.38 days. Third, although the MHD model is 
time-dependent, we assume that the flow at the inner boundary is time-stationary. 
Thus the flow is “corotating", so that rotating spatial variations are responsible for the 
generation of dynamic phenomena in the solution. 

We illustrate this technique in Figure 1 for the time period August 10 to September 
8. 1996. This interval coincided with the campaign known as “Whole Sun Month. 
(WSM) occurring just four months after the termination (minimum) of solar cycle 
22. and has been the focus of considerable research (see papers in the special issue of 
Journal of Geophysical Research. May. 1999). In panel (a) we have traced coronal field 
lines to deduce whether they are open or closed. Regions in which the field lines are 
closed are shaded gray, and open field line regions are shaded black. In panel (b). we 
have shaded the photosphere according to the algorithm described above. Note that 
the central band, indicated with zero speed, does not contribute to the inner boundary 
conditions of the heliospheric model since the field lines there are all dosed. In panel (c) 
we show the speed profile mapped along field lines to 30 R s . The mapping process has 
yielded a relatively complex flow pattern. In particular, the slow-speed band encircling 
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the equator contains significant spatial structure. 

Using this recipe, we derive the flow speed at 30 R t . Bv imposing momentum flux 
balance at the inner radial boundary, we derive the plasma density Theoretically, there 
is little justification for this assumption: however, Ulysses observations suggest that 
momentum flux is roughly conserved [Riley et ai, 1997], To derive the temperature we 
impose thermal pressure balance. This is again, an ad hoc assumption, however, any 
significant pressure gradients would be quickly minimized by the How. To complete the 
necessary inputs, we use the radial component of the magnetic field. B, directly from 
the coronal solution. 

3. Structure of the Heliosphere at Solar Minimum 

To illustrate the structure of the heliosphere during solar minimum conditions, we 
have modeled the WSM I time period (August 10 to September 8. 1990). Figure 2 
summarizes the large-scale features of the heliosphere using the input flow speed shown 
in Figure 1c. The heliospheric current sheet (inferred from the iso-surface B r - 0) is 
displayed out to 5 AU. A meridional slice of the radial velocity is shown at an arbitrary 
longitude. Blue corresponds to slowest speeds (~ 300 km s“‘) and red corresponds 
to fastest speeds (~ 750 km s' 1 ). Superimposed is a selection of interplanetary 
magnetic field lines, as well as the trajectories of the WIND and Ulysses spacecraft 
in a corotating frame. The structure portrayed in Figure 2 fits well with the general 
picture deduced from solar and interplanetary observations during this time period [e.g- 
Riley et ai. 1999: Linker et ai , 1999]. In particular, the streamer belt showed little 
(< 10=) inclination relative to the heliographic equator, but was deformed northward 
at longitudes of ~ 250 - 300= (at the Sun) due to the presence of an active region. 
Note that in spite of the low inclination of the current sheet, by 4 AU it has developed 
considerable structure, including a fold back on itself. 

Figure 3 summarizes the plasma and field parameters as viewed in a meridional 


plane. The left panel shows color contours of radial velocity, the middle panel shows 
contours of number density (scaled by r 2 . to remove the fall off associated with the 
spherical expansion of a constant speed solar wind), and the right panel shows contours 
of gas. or thermal pressure (scaled by r 3 to remove the 1/r fall off associated with the 
near adiabatic (7 - 1.5) expansion of the wind, as well as the 1/r 2 fall-oir associated 
with the spherical expansion). In each panel, the location of the HCS is marked by the 
solid fine. 

Figure 4 provides a complementary view of the same simulation, showing the 
radial velocity, meridional velocity, azimuthal velocity, scaled number density, thermal 
pressure, and magnetic field magnitude, at a heliocentric, distance of 2.3 AU. In each 
panel, the location of the HCS has been overlayed. The trajectory of a spacecraft located 
at a fixed point in space would trace a horizontal line from right to left in this display. 
Consider, for example, a hypothetical spacecraft located at 2.3 AU and 30= above the 
ecliptic plane. From the top-left panel, the spacecraft would sample predominantly fast 
solar wind for most of a rotation, with a relatively slow speed drop at ~ 90= and a 
relatively fast speed increase at ~ 300 3 . By comparison with the bottom 3 panels, we 
would associate the speed increase with the presence of an interaction region (at ~ 300= 
longitude) and the speed decrease with the presence of an expansion wave (at ~ 60° 
longitude). From the top-middle panel, we deduce that the flow deflections through 
the interaction region are first positive and then negative in the meridional plane. In 
our coordinate system, a positive meridional velocity is equatorward. Thus the flow 
is deflected first equatorward and then poleward. These deflections are in the same 
sense as those observed by Ulysses out of the ecliptic plane and are consistent with the 
outward normal of the interaction region being tilted toward the heliographic equator 
(e.g.. Rizzo and Gosling. 1994]. Similarly, we deduce from the positive, then negative 
flow deflections across the interaction region in the azimuthal plane (top right panel) 
that the outward normal to the interaction region is tilted toward the west (i.e.. into the 
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direction of planetary motion), or. alternatively, the interaction region is approximately 
aligned with the Parker spiral direction 

In Figure 5 we compare 1-hr averaged observations of speed, density, and 
temperature by the Ulysses spacecraft (squares) with simulation results (solid curves) 
during the YVSM interval. At the time, Ulysses was returning to lower latitudes and 
was located at ~ 28° North heliographic latitude, at a distance of ~ 4.3 AU. and on 
the opposite side of the Sun from Earth. Thus the measurements of the line-of-sight 
photospheric magnetic field used to drive the solution were always separated by ~ 180= 
from the actual solar wind that Ulysses observed. The profile at Ulysses consisted of a 
simple, single stream pattern. The model appears to reproduce the essential features of 
the large-scale variations, i.e.. the interaction region and the expansion wave: however, 
the interaction region is not as steep as the observations indicate, and the model does 
not reproduce fluctuations on scales less than a few days. 

In Figure G we make an analogous comparison with WIND plasma data. W IND was 
located in the ecliptic plane at 1 AU and observed a more complex pattern of variations 
that are only partially reproduced by the MHD solution. The simulation predicts two 
fast (> 500 km s' 1 ) streams, of which only one is obviously found in the observations 
(at day 255). The second (at day 240) may be related to the fast stream observed on 
day 243. although it is considerably broader. There are other minor perturbations in 
the model: however, none of them can be reliably matched with WIND observations. 

4. Structure of the Heliosphere during the Declining Phase of 
the Solar Cycle 

To illustrate the structure of the heliosphere during more active conditions, we 
have modeled Carrington rotation (CO) 1892 (day 27-55. 1995), which fell on the late 
declining phase of solar cycle 22. In Figure 7 we compare Ulysses plasma and magnetic 


field observations during its so-called “rapid latitude scan” from September. 1994 (when 
it reached its most southern latitude of S80=), to July, 1995 (when it reached its most 
northern latitude of N80=) with simulation results. From top to bottom, the panels 
compare speed, scaled density, scaled temperature, and the radial component of the 
magnetic field. B, . A two-dav box-car running average has been applied to B r to bring 
out large-scale variations. Although several rotations of the Sun occurred during this 
interval, the most- rapid traversal though the equator occurred during ~ 3 Carrington 
rotations (1892-4). During CR 1892. Ulysses moved from -27.0= to -0.9° and from 
1.43 AU to 1.35 AU. The spacecraft crossed into the northern polar coronal hole on day 
87. 1995 [McComas et al. 1998b]. some 32 days after the end of CR 1892. Thus the 
comparison between simulation and observations is expected to be most appropriate 
during Ulysses’ initial immersion into the slow, dense solar wind. Unfortunately, this 
crossing of the polar coronal hole boundary was complicated by the passage of a coronal 
mass ejection that was embedded within the* slow-fast boundary [Gosling et al. 1995c]. 
The CME structure is not included in the present simulations, which as noted earlier, 
assume that the flow and field pattern are time-stationary in a corotating frame of 
reference. 

In spite of these shortcomings, the overall pattern of interaction regions in the 
simulation bears a strong resemblance to the observations. Although the majority of 
the interaction regions and expansion waves are reproduced in the simulations, their 
strengths do not match: densities are too high in the interaction regions, yet not low 
enough in the rarefactions: speed peaks are not high enough: and temperatures do not 
fall enough. Nevertheless, it. is quite remarkable that the final interaction region and 
northern polar-hole crossing are well reproduced almost two rotations after CR 1892. 
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5. Structure of the Heliosphere Approaching Solar Maximum 

As the Sun moves from solar minimum to solar maximum, the assumption that 
the How pattern at the Sun remains constant for an entire solar rotation begins to 
break down. CMEs and other transient phenomena play an ever increasing role in 
the dynamical evolution of the heliosphere. Nevertheless, we can still speak of a 
“rotationally-averaged” picture of the heliosphere during solar maximum conditions. 
Thus we restrict our usp of solar maximum simulations to describe generic features of 
the heliosphere, and use the models to predict what Ulysses may see when it traverses 
(lie southern and northern poles in 2000 and 2001, iespectiielv. 

Figure 8 shows how the assumed radial velocity profile in the photosphere maps 
into the input conditions for the heliospheric solution. The time period modeled (days 
GG-94 of 1909) corresponds to CR 1947. Figure 8a shows the velocity profile in the 
photosphere, where again, the large expanse of zero speed (purple) corresponds to closed 
magnetic field lines and thus does not contribute to the coronal input conditions. Figure 
81, shows the input speed profile at 30 R s . Comparison of Figure 8 with Figure lb and 
le highlights some striking differences between solar wind flow near the Sun at the 
minimum and maximum of the activity cycle. Whereas at solar minimum the slow wind 
is confined to a relatively narrow band organized about the equator, at solar maximum 
slow solar wind dominates the flow at essentially all heliographic latitudes. This is 
in agreement with the conclusions reached by Wang and Sheelcy 11997). A second 
noteworthy point concerns the way features are mapped out from the photosphere to the 
upper corona. At solar minimum (Figure 1). there is a fairly straightforward mapping 
between the two regions, indicating that the topology of the field lines is relatively 
simple. In contrast, at solar maximum, the relationship appears to be more complex. 
For example, all of the the source of the fast solar wind eminating from the south pole 
at 30 R s maps to a photospheric source that is limited in longitudinal exent. Based on 
the photospheric sources of these flows, it is dear that it would be difficult to predict the 
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morphology of the flows in the upper corona without a model such as the one described 
here. 

Figure 9 presents meridional slices of radial velocity, scaled number density, and 
scaled thermal pressure for CR 1947. The location of the HCS is again superimposed 
(solid line). These images may be compared with those in Figure 3. Several differences 
are apparent. First, the symmetric picture characteristic of solar minimum is lost at 
solar maximum. In the northern hemisphere, for example, we see the interaction of an 
isolated high -speed stream with slower solar wind ahead of it: an interaction region 
develops at the interface, bounded by a forward (ahead) and reverse (behind) wave. 

At lower heliographic latitudes, we see the effects of the interaction of medium-speed 
flow with slower flow. In particular, at - 530° the two strongest mid-latitude streams 
(located at ~ 0 - GO* longitude in Figure 8) are driving an interaction region whose 
outward normal is essentially radial. Finally, the interface between slow and fast solar 
wind at the southern coronal hole drives an interaction region (just visible at 5 AU and 
~ 560") whose outward normal is tilted toward the heliographic equator. 

Figure 10 shows spherical slices of the three velocity componets (top panels) and 
density, thermal pressure, and magnetic field magnitude (bottom panels) for CR 194.. 
These images may be compared with those in Figure 4. Note that the heliospheric 
current sheet extends up to ~ 55“ latitude in the northern hemisphere and ~ 40“ in 
the southern hemisphere. As we have noted earlier, the pattern near solar maximum 
is much more complex. At. and near, solar minimum, the interaction and rarefaction 
regions develop a pattern of opposing tilts in each hemisphere [e.g.. Pizza and Gosling. 
1994], so that the fronts run from low latitudes in the east to high latitudes in the west 
in both hemispheres. In contrast, at solar maximum, any orientation of interaction 
regions are possible, as is particularly evident from the color contours of number density 
(bottom- left panel) and thermal pressure (bottom -middle panel). Note that the 
small-scale peaks in density and pressure in Figures 4 and 10 arise from the limited 
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spatial resolution in the calculation. A limited number of higher-resolution simulations 
were run to demonstrate that: (1) such features disappear when sufficient resolution is 
present: and (2) the results are not qualitatively affected by the limited resolution. 

6. The HCS at Different Phases of the Solar Cycle 

The HCS is a convenient feature to illustrate the variability of the heliosphere 
during the solar cycle: being passive, it acts as a tracer for the macroscopic structure 
of the heliosphere [e.g., Riley ft ai. 1990], In Figure 11. we use the simulation results 
described here to illustrate the evolution of the HCS during the course of a solar cycle. 
The first panel shows the HCS near solar minimum The single fold in the surface maps 
back to ~ 180 - 270° longitude. From Figure 1. it can be seen to correspond to the 
northward extension of the slow flow band. Note that this fold is radially asymmetric, 
with its outer edge being sharper than its inner edge. Refering to Figure 1. fast solar 
wind at ~ 280“ overtakes slower wind ahead accelerating it and steepening the HCS 
profile. Conversely, at the inner edge of the deformation (~ 180’ longitude in Figure 1). 
the HCS is “stretched out” as faster wind outruns slower wind creating a rarefaction 
region. It is likely that this picture will differ from the ones produced by kinematic 
models of the HCS. which assume that the speed of the plasma remains constant as 
the plasma moves out from the Sun [e.g.. Sanderson et ai. 1999]. In the middle panel, 
the HCS during the time of Ulysses’ “rapid latitude scan" is shown. Now two folds 
dominate the HCS profile and are in qualitative agreement with the picture built up by 
Smith et al. [1995] based on in situ crossings of the HCS by Ulysses. Finally, we present 
the HCS for CR 1947 in the bottom panel. The HCS is again dominated by two folds 
that extends up to ~ 40 - 55°. In contrast to solar minimum, however, the folds are 
more symmetric, since the wind surrounding it does not have the same large variations 
in speed. 


7. Summary and Discussion 

Global MHD models provide a useful contextual basis with which to interpret 
coronal and solar wind observations. In this study, we have introduced a 3-D resistive 
MHD model of the solar corona and heliosphere, driven by the observed line-of-sight 
photospheric magnetic field. We have illustrated its use by simulating three specific 
time periods corresponding to solar minimum, the declining phase of the solar cycle, 
and near-solar maximum conditions. We find that the model reproduces the essential 
features of the observations at solar minimum and during the declining phase. 

The near-solar maximum simulations suggest lower, more variable speeds and the 
absence of fast streams over a much broader range of heliolatitudes, in agreement with IPS 
measurements and recent Ulysses observations (http://swoops.lanl.gov/recentvu.html). 
Our simulations suggest that the pattern of alternating tilts that dominated Ulysses 
observations during its first orbit will likely be lost at solar maximum. To produce this 
pattern of “alternating tilts”, slow solar wind must emanate from lower latitudes than 
fast wind. Near solar minimum, this is accomplished because a band of slow solar wind 
is organized about the solar equator. Either a simple tilted-dipole band [Rizzo. 1991] 
(where the wind is organized about the magnetic equator which is tilted with respect 
to the rotation axis), or a warped band (where a warped flow is organized about the 
rotation axis) is sufficient to produce a regular pattern of alternating tilts. On the other 
hand, when isolated coronal holes (i.e.. regions of faster wind completely encircled by 
slower wind) are present, the simple tilted fronts are replaced by more complex shapes. 
In the simplest case, a spherical coronal hole at low- or mid-latitudes would produce 
a U-shaped interaction region in the solar wind. Depending on a spacecraft’s latitude 
relative to the center of the disturbance, it would see a variety of orientations. 

Based on the results presented here (and others not shown), we can make some 
predictions as to the nature of the solar wind that Ulysses will encounter during its 
crossings of the polar regions. First, the simulations suggest that while Ulysses may 
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encounter fast (i.e.. ~ 750 km s ') streams, it will likely not remain in them for more 
than a fraction of a solar rotation. This is in agreement with the predictions made bv 
HTim? and She.elcy [1997] . As discussed above, although the pattern of alternating tilts 
at Ulysses will likely disappear at low- and micl-heliographic latitudes, if a polar coronal 
hole is present, this pattern may return. Finally, the results suggest that crossings of 
the HCS may persist up to ~ 70° heliographic latitude. 

Acknowledgments. The authors gratefully acknowledge the support of the National 
Aeronautics and Space Administration (Sun- Earth Connections Theory Program. Supporting 
Research and Technology Programs, and SOHO Guest Investigator Program) and the National 
Science Foundation (Space Weather Program ATM9613834) in undertaking this study. We 
also thank National Science Foundation at the San Diego Supercomputer Center and the 
National Energy Research Supercomputer Center for providing computational support. 


References 

Axford. W. I., Interaction of the solar wind with the interstellar medium, in Solar Wind, 
edited by C. P. Sonett, P. J. Coleman, and J. M. Wilcox, NASA Spec. Pubi. NASA 
SP- 308. 609. 1972. 

Bame, S. J., D. J. McComas, B. L. Barrarlough. .1. L. Phillips, K. J. Sofaly. J. C. Chavez, B. 

E. Goldstein, and R. K. Sakurai. The Ulysses solar wind plasma experiment, Astron. 
and Astrophys.. Suppl- Scr., 92. 237. 1992. 

Balogh, A.. T. .1 Book, R. J. Forsyth, P. C. Hedgocock, R. J- Marquedant, E. J. Smith. D. 

J. Southward, and B. T. Tsurutani, The magnetic field investigation on the Ulysses 
mission: Instrumentation and preliminary scientific results. Aslrcm. Astrophys, 92, 221, 
1992. 

Feldman, W. C„ B. L. Barradough. J. T. Gosling. D. J. McComas. P. Riley. B. E. Goldstein, 
and A. Balogh, Ion energy equation for the high-speed solar wind: Ulysses observations. 

J. Geophys. lies.. 109, 14,547, 1998. 

Fisk. L. A.. Motion of the footpoints of heliospheric magnetic field lines at the Sun: 

Implications for recurrent energetic particle events at high heliographic latitudes, J. 
Geophys. Res.. 101, 15,547, 1996. 

Gosling. J. T„ S. J. Bame. D. J. McComas. J. L. Phillips. V. J- Pizzo, B. E. Goldstein, and M. 
Neugebaucr, Latitudinal variation of solar wind corotating stream interaction regions: 
Ulysses, Geophys. Res. Lett.. 20, 2789. 1993. 

Gosling. J. T„ S. J. Bame, D. J. McComas, J. L. Phillips. V. J. Pizzo, B. E. Goldstein, and M. 
Neugebaucr, Solar wind corotating stream interaction regions out of the ecliptic plane: 
Ulysses. Space Sci. Rev., 72. 99. 1995a. 

Gosling. J. T.. W. C- Feldman. D. J. McComas, 3.L. Phillips. V. J. Pizzo. and R. J. Forsyth. 
Ulysses observations of opposed tilts of solar wind corotating interaction regions in the 
northern and southern solar hemispheres, Geophys. Res. Utt., 22. 3333. 1995b. 

Gosling. J. T.. S. J. Bame, W. C. Feldman, D. J. McComas, .1. L. Phillips, B. E. Goldstein, 
M. Neugebaucr. 3. Burkepile, A. J. Hundhausen. and L. Acton. The band of solar wind 



21 


22 


variability at low heliographs latitudes near solar activity minimum: Plasma results 
from the Ulysses rapid latitude scan. Geophys. Res. Lett., 22, 3329, 1995c. 

Gosling, J. T.. S. J. Bame. W. C. Feldman. D. J. McComas. P. Riley, B. E. Goldstein, and M 
Ncugebauer. The northen edge of the band of solar wind variability: Ulysses at ~ 4.5 
AU. Geophys. Res. Lett.. 24. 309, 1997. 

Jokipu. R. J.. and B. Thomas, Effects of drift on the transport of cosmic rays. IV. Modulation 
by a wavy interplanetary current sheet. Astrophys. J.. 243, 1115, 1981. 

Levine. R. J.. M. D. Altschuler, and J. W. Harvey, Solar sources of the interplanetary magnetic 
field and solar wind. .7. Geophys. Res., 82. 1061, 1977. 

Linker. .1. A., and Z Mikic. Extending coronal models to Earth orbit, in Coronal mass 

ejections. Geophys. Monograph 99. edited by N. Crooker. J. Joselyn. and J. Feyntnann, 
American Geophysical Union, Washington, D.C., 269. 1997. 

Linker. .1. A.. Z. Mikic. D. A. Bisecker, R. J. Forsyth. S. E. Gibson, A. J. Lazarus, A. Lecinski. 
P. Riley. A. Szabo. and B. J. Thompson, Magnetohydrodynamic modeling of the solar 
corona during whole sun month. J. Geophys. Res. 104. 9809, 1999. 

Lionello. R., Z. Mikic. and J. A. Linker. Magnetohydrodynamics of solar coronal plasmas in 
cylindrical geometry. J. Comput. Phys.. 140. 172. 1998. 

Lionello. R.. .J. A. Linker, and Z. Mikic. Including the transition region in models of the 
large-scale solar rorona, to appear in A s trophy s. J., January. 2001. 

D. J. McComas, S. J. Bainc, B. L. Barraclough, W. C. Feldman, H. O. Funstcn. J. T. Gosling. 
P. Riley. R. Skoug. A. Baiogh. R. J. Forsyth, B. E. Goldstein, and M. Ncugebauer, 
Ulysses return to the slow solar wind. Geophys. Res. Lett., 25. 1. 1998a. 

D. J. McComas, P. Riley. ,1. T. Gosling, A. Baiogh. and R. J. Forsyth, Ulysses' rapid crossing 
of the polar coronal hole boundary. J. Geophys. Res.. 103. 1955, 1998b. 

Mikic, Z., and J. A. Linker, Disruption of coronal magnetic field arcades. Astrophys. J., 430. 
898, 1994. 

Mikic. Z., and Linker. J. A., The large-scale structure of the solar corona and inner heliosphere, 
in Solar Wind Eight, edited by D. Winterhalter. J. T. Gosling, S. R. Habbal. W. S. 
Kurth, and M. Ncugebauer. published by American Institute of Physics, 1996. 


Mikic, Z„ J. A. Linker, D. D. Schnack. R. Lionello. and A. Tarditi, Maguetohydrodynaniic 
modeling of the global solar corona. Phys. Plasmas, 6. 2217, 1999. 

Ogilvie. K. W.. et ai., A comprehensive plasma instrument for the WIND spacecraft. Space 
Set. Rev.. 71. 55. 1995. 

Pizzo. V. J.. Interplanetary shocks on the large scale: A retrospective on the last decade's 
theoretical efforts, pp. 51-68, Geophys. Monograph 35. published by the American 
Geophysical Union, Washington, D. C., 1985. 

Pizzo. V. ,1., The evolution of corotating stream fronts near the ecliptic plane in the inner solar 
system, 2, Three-dimensional tilted-dipole fronts, J. Geophys. Res.. 96. 54(15, 1991. 

Pizzo. V. J.. Global, quasi-steady dynamics of the distant solar wind, 1. Origin of north-south 
flows in the outer heliosphere. J. Geophys. Res., 99, 4173. 1994a. 

Pizzo, V. J.. Global, quasi-steady dynamics of the distant solar wind. 2, Deformation of the 
heliospheric current sheet. J. Geophys Res.. 99. 4185. 1994b. 

Pizzo. V. .1., and J. T. Gosling, 3-D simulation of high-latitude interaction regions: Comparison 
with Ulysses results. Geophys. Res. Lett.. 21. 2063, 1994, 

Riley, P„ J. T. Gosling, L. A. Weiss, and V. J. Pizzo, The tilts of corotating interaction regions 
at mid-heliographic latitudes, J. Geophys. Res.. 101. 24.349, 199G. 

Sanderson, T. R.. D. Lario, M. Maksimovic, R. G. Marsdcn, C. Tranquille, A. Baiogh. R. J. 
Forsyth, and B. E. Goldstein. Current sheet control of recurrent particle increase at 4-5 
AU. Geophys. Res. Lett., 26. 1785, 1999. 

Sarabhai. V., Some consequences of nonuriifonnity of solar wind velocity. J. Geophys. Res.. 

68. 1555. 1963. 

Smith, E. J.. and J. H. Wolfe, Observations of interaction regions and corotating shocks 
between one and five AU: Pioneers 10 and 11, J. Geophys. Res.. 3. 137, 1976. 

Smith. E. .1.. A. Baiogh. M. E. Burton. G. Erdos, and R. J. Forsyth. Results of tl.e Ulysses 
fast latitude scan: Magnetic field observations, Geophys. Res. Lett.. 22. 3325. 1995 

Totten, T. L., .]. W. Freeman, and S. Arya, Application of the empirically derived polytropic 
index for the solar-wind to models of solar- wind propagation, ./. Geophys. Res.. 101. 
15,629. 1996. 



23 


21 


Tsurutani. B. T.. and E. J. Smith, Interplanetary discontinuities: Temporal variations and 
radial gradient, J ■ Geophys. Res.. 84. 2773. 1979. 

Wang. Y. -M.. and N. R. Slioeley. Jr.. Solar wind speed and coronal flux-tube expansion. 

Astrophys. J.. 855. 726, 1990. 

Wang. Y. M. and N. R. Shveley. Jr.. The high latitude solar wind near sunspot maximum. 

Geophys. Res. Lett.. 24- 3141. 1997. 

"TX Linker. Z. Mikic. and P. Riley. Science Applications International Corporation. 

San Diego. CA 92121. (e-mail.linktTairis023.saic.rom: iiiikicUiris023.saic.com.uk2yhaveti.saic.com) 

Received N/A: revised N/A: accepted N/A. 


Figure 1. Figure 1; Derivation of the ratlial speed profile at 30/?,: (a) Location of 
open and closed magnetic field lines in the photosphere: (b) Radial speed profile based 
on equation (7): and (c) Speed profile at 30/?,. 

Figure 2. Model solution for CR 1912-1913. The heliospheric current sheet (inferred 
from the iso-surface B, = 0) is displayed out to 5 ALL The central sphere marks the 
inner boundary at 30/?,. A meridional slice of the radial velocity is shown at an arbitrary 
longitude. Blue corresponds to slowest speeds (~ 350 km s' 1 ) and red corresponds to 
fasted speeds (~ 750 km s ') Superimposed is a selection of interplanetary magnetic 
field lines originating from different latitudes. Finally, the trajectories of the WIND and 
Ulysses spacecraft are marked, 

Figure 3. Large-scale structure of the heliosphere in the meridional plane for the same 
interval as Figure 2. Left panel shows color contours of radial velocity, middle panel 
shows contours of number density (scaled by r 2 ). and right panel shows contours of gas. 
or thermal pressure (scaled by r 5 ) at an arbitrary longitude. In each panel, the location 
of the IK'S is marked by the white curve. 

Figure 4. From top left, in a clockwise direction: the radial velocity, meridional ve- 
locity. azimuthal velocity, scaled number density, thermal pressure, and magnetic field 
magnitude, at a heliocentric distance of 2.3 AU. In each panel, a black line maiks the 
location of the HCS. 

Figure 5. Comparison of Ulysses in situ measurements of speed, density, and tempera- 
ture (squares) with simulation results (solid curves) for WSM time period. 

Figure 6. Comparison of WIND in situ measurements of speed, density, and temperature 
(squares) with simulation results (solid curves) for WSM time period. 

Figure 7. Comparison of Ulysses observations (blue) and simulation results (red) of 
speed, density, and temperature from September. 1994 to July. 1995. 

Figure 8. (a) Speed profile in the photosphere for time period day GG-94. 1999. corre- 
sponding to CR 1947. (b) Inferred speed profile at 30 It,. 
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Figure 9. Meridional slices (at an arbitrary longitude) of radial velocity, scaled number 
density, and scaled thermal pressure for CR 1947. The white curves mark the location 
of the HCS. 

Figure 10. Spherical slices of the three velocity componets (top panels) and density, 
thermal pressure, and magnetic field magnitude (bottom panels) for CR 1947. 

Figure 11. HCS at the three epochs of the solar cycle described in the text: (a) solar 
minimum: (b) the declining phase: and (c) approaching solar maximum. 
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Abstract. We present a 2.5D axisyminetric MHD model to self-consistent ly 
describe the formation of a stable prominence that supports cool, dense material in 
the lower corona. The upper chromosphere and transition region are included in the 
calculation. Reducing the magnetic flux along the neutral line of a sheared coronal 
arcade forms a magnetic field configuration with a flux-rope topology. The prominence 
forms when dense chromospheric material is brought up and condenses in the corona. 
The prominence sits at the base of a helmet streamer structure. The dense mateiial is 
supported against gravity in the dips of the magnetic field lines in the flux rope. Further 
reduction in magnetic flux leads to an eruption of the prominence, ejecting material into 
tiie solar wind. 
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1. Introduction 

Prominences (called filaments when observed on the solar disk) support cool, dense 
chromospheric material (~ KPK and 10 ,u - 10“ /nil 3 ) against the solar gravity in the 
surrounding hot. tenuous corona (~ 10 6 K and 10 7 — It) -1 /cm 3 ). They are observed to lie 
above magnet ie neutral lines in the photosphere and near the base of helmet streamers 
(regions of closcxl magnetic field that have confined the coronal plasma). The magnetic 
field in the prominence often exhibits ‘inverse polarity.” meaning that when the c oronal 
magnetic fields embedded in the prominence cross over the neutral line they point in 
the direction opposite to that indicated by the photospherie magnetic field polarity. 

The prominence magnetic field is itself nearly aligned with the filament channel [ Martin 
ct ai. 1994: Martin and Echols. 1994]. indicating a highly sheared (and therefore 
magnetically energized) configuration. 

Prominences haw been studied for many years, yet the means by which these 
structures form and are maintained is still not understood, nor is their violent eiuption. 
Three main difficulties confront any prospective theory attempting to describe the 
formation and evolution of prominence’s: (1) Finding a magnetic configuration with 
“dips' (concave upward portions of fluxtubes) that can gravitationally support the 
dense material: (2) understanding the mechanism hv which chromospheric material is 
trapped in the clipped field lines and maintained there to form a condensation: (3) 
elucidating the process that leads to the release of magnetic energy and the disruption 
of these structures. Bec ause of the great complexity of the entire problem, these’ three 
individual aspects have typically been approached separately. 

Models of magnetic: field configurations for prominence support that develop 
the required dips and inverse polarity usually compute force-free magnetic fields 
and assume that the prominence material provides only a small perturbation to the 
magnetic structure [eg., ran Dalkgooijen and Martens. 1989. 1990: Antiochos ct ai. 1994; 
Aulanicr and Dcmonlin. 1998: ,4 man <1. al. 1999] Given a favorable structure for 


supporting the filament, computation of the complicated dynamics and thermodynamics 
of condensations have been performed by assuming that the plasma flows along fixed 
magnetic flux tubes, which reduces the problem to one-dimensional hydrodynamics with 
energy transport [ Poland and Mamsku. 1986; Mok ct ai 1990: Antiochos and Klimchvk 
1991: Antiochos ct. ai. 1999a]. 

Models of prominence eruption typically start from configurations favorable for 
prominence support [e g., van Dallcgooijen and Martens. 1989: Priest and Forbes. 1990; 

I sen berg ct ai. 1993) and are closely related to the problem of coronal mass ejection 
(CME) initiation [Forks and Priest. 1995: Linker and Mikic. 1995: Low. 1997; Miku: 
and Linker. 1997: IPn and Guo. 1997; Antiochos el ai. 1999b: Lin and Forbes. 2000], as 
these phenomena are linked obsei vat ionally \Hundhauscn. 1997] and require the release 
of stored magnetic field energy and the opening of previously closed magnetic field 
regions [eg. Aly. 1984; S tarrock. 1991: Forbes . 1992: Mikic and Linker. 1994; Antiochos 
1998). In recent simulations [Mikic ct ai 1999: Aman ct ai 1999, 2000] we have found 
that magnetic flux emergence and cancellation in the photosphere can lead to the 
formation of magnetic flux ropes in sheared or twisted arcade configurations. When 
the flux cancellation readies a critical threshold, the entire configuration erupts with a 
considerable release’ of magnetic energy. 

The theoretical investigations described above decouple the salient processes and 
focus on modeling individual aspects of the problem. This approach is useful foi 
revealing the basic underlying physics. However, a complete pic ture of prominence 
formation, evolution, and eruption ultimately requires a comprehensive model of all of 
the processes together. This is partic ularly true now that different models can more or 
less equally describe the basic features of the observations, albeit somewhat superficially 
Eventually, models will need to produce more detailed predictions that can* be tested 
direc tly by observations (for example, by producing simulated emission that can be 
compared with images from spacecraft such as SOHO or TRACE). This challenging 



goal requires that the complex thermodynamic processes of the ti])per chromosphere and 
transition region lie incorporated into multi-dimensional magnetohvdrodyuaniic (MHD) 
computations. 

The purpose of this paper is to demonstrate that we can now begin to study the 
problem of prominence formation and eruption with a more comprehensive approach 
We show that when energy transport processes are included into calculations similar 
to our recent models of prominence support and eruption, chromospheric- like material 
can be trapped on helical field lines and lifted into a stable' configuration in the lower 
corona as a result of flux cancellation at the neutral line. When further flux cancellation 
occurs, the entire configuration erupts into the corona. 

The plan for the rest of the paper is as follows: In the next section (2). we describe 
our computational mcthexlology. Section 3 describes results for both poly tropic and full 
thermodynamic simulations, and section 4 discusses the implications of our work and 
directions for future progress. 

2. Methodology 

Lionello et at. (1999a: 2001] describe the use of our time-dependent 3D MHD model 
to compute helmet streamer configurations in the solar corona that include the lower 
solar atmosphere (upper chromosphere and transition region) Those calculations form 
the basis for the studies discussed in this paper: here we note some important details 
relevant to the present calculations. 

In this paper we confine ourselves to azimuthallv symmetric, two-dimensional 
solutions (with three components of magnetic field and velocity). Wo solve the following 
set of equations in spherical coordinates: 



r 


VxD 

1 dB 

rl)i 


-V x E. 


( 2 ) 


fi 


E -I -v x B =■ i/J, (3) 

r 

^ + V-(pv) - 0. (4) 

()i 

~“T {iJT + V ' VT ) “ ~ 7 ' V ’ V - ]~( V ' + "rTI„Q(T) - //„. - tl„ - //„). (5) 

/> + v-Vv) — -J x B - Vp + pg + V • (izpVv). (fc>) 

where B is the magnetic field. J is the electric current density. E is the electric field. 
j>. v. 1 1 . and T are the plasma mass density, velocity, pressure, and temperature, 
g - —gufli^/r 2 is the gravitational acceleration (with /f, the solar radius), i] the 
resistivity, and v is the kinematic viscosity. In the energy equation. Eq (5). Q(T) is the 
(optically thin) radiation loss func tion as in Athay [1980]. n, and n r are the electron 
and proton number density (which are equal for a hydrogen plasma). ~ * 5/3 is the 
poly tropic index. // rh is the coronal heating term (a parameterized function). If ri = y.l 1 
and //„ = ft tzVv : Vv are the ohmic and viscous heating terms (neglected in these 
simulations), and q is the heat flux 

Coronal heating is specified as a nommiform profile that varies exponentially 
with radial distance, and both the' heat deposition length scale (A) and flux vary with 
latitude The length scale varies from A - 0.7 fis at the poles to A -- 0.1 Its at the 
equator, and the heat flux at r = lis is 10 s erg cm 2 s~ 1 at the poles and 5 x 10 1 erg 
cm 2 s 1 at the equator. These parameters have been chosen to yield the strongly 
concentrated heating and higher densities seen in active regions near the equator of the' 
simulation while also providing the more distributed heating present in coronal holes in 
the open field portions of the simulation. Liondlo f t at. [200 1 j describes further details 
of the coronal heating parameters and the use of both eollisional (Spitzer s law) and 
collisionless ( Holtwrg. 1978] thermal conduction as a function of radial distance from the 
Sun. For the results shown in this paper, we have modified the Spitzer law to reduce the 
steepness of the temperature and density gradients in the lowest part of the transition 



7 

region (see the Appendix for a more detailed discussion). We have found that this 
procedure allows coarser ineshes to be used and yields solutions that are qualitatively 
similar to those obtained with the full Spitzer thermal conductivity. 

The calculation described here has been performed on a 201 x 201 nonuni form 
(r. 0) grid, with the mesh points highly concentrated in the equatorial region near the 
lower boundary. A mesh with A r « 6 x 10 ‘'R, in the upper chromosphere (near the 
lower radial boundary) and A 0 k. 0.2° near the equator is used. Other calculations 
with 121 x 101 mesh points, with coarser mesh resolution, were also performed prior 
to this calculation in order to scope out the parameters and find a regime in which a 
prominence was formed. The results with the higher- resolution mesh allowed us to use 
a larger Lundquist number, and we report only these results here. This simulation used 
20 hours of CPU time on a Cray TOO supercomputer, and we did not deem it necessary 
to pursue calculations with higher mesh resolution considering the limitations on our 
computing resources. We believe that increasing the resolution beyond that used in the 
present calculation would not change the' results significantly. However, when we begin 
modeling more realistic prominences in the future we will verify that a further increase 
in the spatial resolution does not change the results. 

The simulation domain extends out to 30/V,,. At this upper boundary the flow 
is supersonic and supor-Alfvenic-. and we implement boundary conditions that utilize 
the characteristics to allow only outgoing waves there [ Linker and Mikic. 1997], 

A uniform resistivity V has been used, corresponding to a resistive diffusion time 
th — 47r/?, 2 /(r/f J ) = 4 x 10 1 hours (for a length scale of /?,). Underneath the relaxed 
helmet streamer, the Alfven speed (l',t 0 ) near the equator at 10.000 - 20. 000 km altitude 
is about 1830 km/s. so the Alfven travel time (r/i = R»/V. w) is 0.3 minutes, and the 1 
Lundquist number Tu/t,.\ = 3 8 x IIP. A uniform viscosity u is also used, corresponding 
to a viscous diffusion time r„ = R» fv such that t„/t a = 3.8 x 10 s . Typically, we 
find that much higher values of r« can be specified for our algorithm if t„ is kept at a 
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relatively smaller value. 

The method of solution of (1-0) lias been described previously [Mikic and Linker. 
1994: Linker und Mikic. 1997: Lionclln e.t al. 1999b: Linker ct at. 1999: Mikic f t 
al. 1999]. For the parameter regime of the calculations wo describe here, stiffness of the 
equations is introduced by the combination of both the high Alfven speed and thermal 
conduction in the transition region and lower corona, and the use of small niesli cells 
in this region to capture the steep gradients in density and temperature there. The 
stiffness introduced in the time-integration of the equations by the Alfven speed is 
treated efficiently using a serni-implieit method [Mikic and Linker. 1994). The accurac y 
of the semi-implicit method has been studied previously [Schvack rt al.. 1987: Mikic 
ct al.. 1988], and it is known that when the time step falls below the CFL limit, the 
semi-implicit method is identical to an explicit method. The stiffness introduced by the 
parabolic equations resulting from the thermal conduction, resistivity, and viscosity, are 
treated using standard fully-implicit methods. 

3. Results 

Aman ct al. [1999; 2000] have described how magnetic flux changes at the 
photosphere, when introduced into a three-dimensional sheared arcade field, c an yield 
solutions with stable* magnetic flux ropes suitable for prominence support, and how 
further magnetic flux changes can lead to magnetic energy release and eruption These 
calculations were performed in the "zero beta” limit of equations (1-0): equations (1-3) 
and (0) are solved with P — 0 and a fixed profile assumed for p (equilibria found in 
this manner are force-free solutions). Mikic and Linker (1999] have shown that when 
the same procedure is applied to polytropic MHD solutions with helmet streamer 
configurations [solutions of equations (1-6) with the energy source terms in (5) set to 
zero and 7 — 1.05]. the streamer disrupts and material is ejected out into the solar wind. 
Figure 1 shows the eruption of a flux rope in a three-dimensional calculation of this kind. 



A detailed description of these results in presently in preparation. Here we describe 
how the inclusion of processes in the lower solar atmosphere in these calculations [full 
solutions of equations (l-6)j leads to the* formation of a prominence-like structure and 
its subsequent eruption 

3.1. Helmet Streamer Solution 

In the first phase of the* calculation, we generate a helmet streamer equilibrium [e g.. 
Linker and Mtku':. 1995]. We start with a potential magnetic held in the corona that 
matches a specified distribution of radial magnetic- field at the solar surface li rV . The 
distribution we choose is the sum of a weak dipole field (/* r0dlp = 2.2G af the poles), 
and a stronger {H rilbip = 9.5G) conc entrated bipole near the neutral line. For simplic ity 
we choose a configuration that is symmetric in latitude* about the* solar equator. Our 
intent in including the strong hipole is to model the large-scale effect of the field in 
an active region (where* prominences frequently form). The equatorial position of the 
bipole at the equator is not very realistic for prominences, but is convenient for the 
illustrative purposes of this calculation. We place our lower boundary at the top of the 
chromosphere, and we impose a fixed temperature T 0 = 20.000K and a plasma number 
density n — n p = n, — llfl'em -3 . We impost* a spheric-ally symmetric solar wind solution 
which includes the* upper chromosphere and the transition region and is consistent 
with the chosen houndary values for temperature and density This combination is not 
initially at equilibrium. We integrate the time-dependent MHD equations in time until 
the solution settles down to an equilibrium (for 3807, , where t* is the* Alfven time 
described in section 2). The final state* models the coronal plasma, and includes the 
transition region in the* calculation. This procedure hies been described by Linnrllo it 
(il. [2001], The final solution has a coronal streamer with c losed field lines, surrounded 
by open field line's along which the solar wind Hows outward. Figure 1 of Linker and 
Mik/e [1995] shows an example of the magnetic configuration of a 2D helmet streamer 
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using a polytropic energy equation; Figure 1 and Figure 3 of Linnrllo cl ai [2001] show 
example's of a helmet streamer configuration computed with thermodynamics. 

3.2. An Energized Helmet Streamer 

In the second phase* of the* simulation, we* apply a shear flow near the neutral line* 
that builds free magnetic energy into the streamer. This shear flow is not intended to 
model ac tual flows on the Sun. It is just a convenient mechanism for producing strongly 
sheared field lines that are nearly aligned with the neutral line, a frequently observed 
characteristic of filaments [Martin and Echols. 1994). We use a shear profile* that is like 
the one used by Miktc and Linker [1994] . with a width A () m = 8.5°. The shear is applied 
from t = 380t,., to t = 570 t a . with a maximum velocity r 0 = O.OOfiEw Figure 2a and 3a 
show contours of the flux function (projec tions of the magnetic field lines) superimposed 
on the plasma density (Figure 2a) and plasma temperature (Figure 3a). at the end of 
tlu* shearing phase. The maximum displacement of the footpoints from their original 
position is 0.8 R„. and the sheared field has a magnetic energy equal to 2.3H' P „ ( . The 
calculation is them continued in a relaxation phase to ! ^ 589r,,. 

3.3. Creation of a Magnetic Flux Rope and Subsequent Eruption 

lu the final phase of the calculation, we change the magnetic flux at the boundary 
to generate* a flux rope. On the Sun it is frequently observed that the magnetic fields 
in an active region tend to disperse days to weeks after its emergence. During this 
time, filaments are frequently observed to form along the* neutral line. At times, these* 
filaments disappear, presumably due to eruption, and may even reform in the* same 
location later. T his dispersal of magnetic flux is thought to occur on a small spatial scale 
by annihilation and submergence of magnetic dipole elements, and has been modeled as 
a convective-diffusive process on a large scale [Wang ct ai. 1989: Wang and Shcc.lcy. 
1990] The disappearance of photospherie magnetic flux has long been suspected of 
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MX lines of the flux rope The height of the sunulu.nl pro.. f» »>■»»' 15 -™ 

k,„ As .Uo calculation is oomim.o.l further II,.' flux ropo slowly .l.lfns.'s booauso of .1,0 

presence of finite resistivity. 

If flux reduel ion is halted at » lovol of 7 r« of .ho initial amount »'“1 ,h «' - ah ula.lon 
is continued, the flux ro,,e and helmet streamer erupt in a manner similar to the ease 
shown in Fiflure M. We have found that the exact level of flux reduction nocossarv for 
eruption in a ttivet, configuration depends on the details of the surface flux rhstr.l,. 
and the amount of shear introduced. Prior to reaching the critical level, the behavior of 

thi' system is quasistatic. 

Figure 5 shows the prominence height as a function of time during the eruptive 
phase. The eruption that was prodneed ill this ease was not very fast (reaching - 40km/s 
a, 211.). although it is still aeeelera.ing at this stage This is in .outlast t. some of the 
poty, topic simulations we have ml, with this mechanism (speeds > a#km/s). A mote 
concentrated shea, and a larger initial magnet, c field strength (enabling release <* "tore 
magnetic energy) may be , cubed to create a fas, mass ejection containing a filament. 

4. Discussion 

The possibility that prominences are supported by magnetic flux ropes has been 
considered in recent year, by a number of authors (e g. nan ****** mt *•**• 
t titty, 199(1. /last and Kumar. 1994, Cl ten. 1996. A, .lamer and Demo, dm. 1998|. On, 

work \Amm H «l- 1™: W») «* "* "**« ” h “ VC l”*"" 1 l '"" 
many respeels to rite force-free calculations of run Balk, m a, and Martens |1989(. 
which showed bow flux cancellation could lead to the formation of flux rope*. The 
eompiltations we have described here show ,l,„ when the full MHD equations are used, 
including an energy equation that takes into account energy transport in the ttppe, 
chromosphere and transition region, chromospheric material can he lilted by the lit lit al 
field Hues and supported against gravity inside a helmet streamer The flux rope can 
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tutlrtrq.tcullv he ejected into tile sola, wind as a consequence of eruption. 

The idealised calculations presented here are intended to illustrate that h is possible 
to perform sclf-cousistent MHD ralculalions of prominences. Clearly, a 2D axisyimuemt 
model has limited applicability to solar observations, hr rhts geometry the helical field 
lines of the flux rope form a toms around the Sun The connection of the magnetic 
fields embedded in the ptommenoe to the photosphere, and the ,««sihte role of siphon 
flow, in hlhns/draimng I he prominence. ,s an important question mb addressed by Has 
calculation. Preliminary onmtUmensional solutions for energy flow along the magnet.e 
field loops from a three-dimensional flux-rope configuration indicate that condensations 
,1„ indeed form, but a fully self-ronsisleni 3D simulation is required in investigate tins 
question more completely. We have po, formed 3D simulations of flu*-ro,.c formation in 
P„tv, topic helmet streamers (Figure 1). and we intend to extend the present results to 

3D as well. 

Once the simulations become more realistic, it will l,e necessary to explain the 
many detailed (ea.nrtw that have been observed in prom.ncnccs (e.g. Martin el * 
1994 | For example, what are the impheations of the olrscrved relationship between the 
filament axtal field and the skew of eorottal arcades [Martin and MeAllmler. 19 * 1 . and 

relationsliip l.e.weot, ddferenttal solar rotation and the axial field in polar crown 
filaments jean Da llw« ei al. 19981? Do these considerations imply that the axtal 
field originates front Mow the photosphere \Prust el al 1996]? Magnetic flux rope, 
are no, the only candidates lot explaining prominence support (Martin and Eclml,. 
,994- Anttoefws e, al 1994) and eruption [dnimeta e, «1. 19991,]. Unrletstandmg the 
pte-eruptive state of filaments, as well as .holt violent e, upturn, tequues ealeulattous 
that ran pmlic, nlscr.ablc qua, flit, e». We consider the calculations «f l*~ I"**" 1 
here as a hist step of that process. 
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Appendix: Modified Thermal Conduction 

,» ,1„. lower mmsili,,,, n*«m. tbm hr . ctar balance '«««• 

(V .q in oi„„,i,m (Ml comms ta» »«- '»*■ <<mma "’" l 

away' in (Mi by -he .««*•» "™»' <”> J "*> 

thr Spit zer form of the thermal eondurtivitv is appropriate: 

q — — K,|bb • VT. (A1) 

- „TV (« = !) x lO 7 ill I S-", units) . b in 111.' .mil vetor along .be W tk 
field, and T is the temperature in K 

Tin- Minimus of romponiture in the n-Misitiim h™ ll “' 

of conductive M flux fro... tin- two* -1 «mtw l—> hr the 

The iomporaturo dependence of Si' 1 '"" ^ * 7 "''’ 

mnperainre grailiii.t in tl.o W » ansi. io.i rogion. Urn*. « aL j200I| used 

„ig„lv resolved mosi,,, (A, . 8 x W*H.. « IK km) .0 emtae «* " ™ 

an nr.ifi. ial ImM, of this gradient should be ta** »»!.<"" ..u.ili.nt.vofy 
, hanging ,1.0 ovorrrll sol,,,,,,!,. With tins goal in rninrl. wo U.vo fi.r,.n,l.,o,l », in sml, » 
way as ,l,o preserve tl.o osseni.al physic of Iho Spiixor model for the npl.or «—>» 
rogion anil corona, wlnlo ...hroing iho siccpnes of the m«K ■» «* kme. !'«' <>' 
transition region: 

«, - «»[xr OT + t)-»)r*JiS"l- lA2) 

*m - ^(lt.anh^^'))- ,AJ ’ 

T,io function . varies smoothly bo, ween 0 mt 1 and rognla.c the ansi., on tan .1,0 

Spi,,,, of «„ (when T i T, ami the mmMed regime («■« T >’ 

,„o m»WW login, o. whioli a, 1 , .lie ,0 the lower iransiiion region. „ X T-. «n,l 


1G 


„ power ,s II < „ < S/2. The iransiiion ho.woon the «. -rf— * 

tempera! me interval ilopcnling on A7;„,.,. For Iho sinmlalion rculls |.ros,-n,o,l lloro. 
wo chore n «= 0 (oi,rrc|Himhng to rolls' ani i M for 1 < Tiinoil). M„. ,1 - 251). 0(H1 K. and 
A7',„„d ~ 20. 000 K. 

Tins «, Wta* »>**" k ’ W '"' 

rogion . allowing t ho, ..ml oonduoii,,,. .0 balance mdtam. k» «■* » 

gradiollt Wo have found .ha. this mod, Hod l<«m of Iho thermal on, nine, ion yields 

solnt ions II, at are .,naht«tivoly similar In those obtained with the full Spinier Hiormal 

i.indiiotivily. and wo wore ahlo ,0 radial grid sis.o in the lower tttmmta. 

rogion l,y almost an or, lor of magniludo over that userl in tamJlo ct of. 12001]. 

Quant itativoly. ,1.0 tempemtnre of solnlions on individual loops maiohos cicely in the 
oornna nllhongh dilfol cnees in the donsilv mil lie as largo as ,1 fa, lor of 
i.V, regard llio moditiod lliormal oonduotion as a nsoinl lis'l for ixploiing tin i| 
properties of solutions without using extremely fine meshes. 
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Figure 1. Simulated eruption of a magnetic flux rope on the Sun. (a) A three- 
dimensional flux rope formed by reducing the- magnetic- flux in a sheared arcade, computed 
with the poly tropic model (see text). The flux rope lic-s within a helmet streamer that 
is surrounded by open magnetic field lines along which the solar wind streams outward 
When the amount of magnetic flux decrease is small enough, the flux rope is stable. 
(h)-(d) Further reduction of the photospheric magnetic field leads to the eruption of the 
flux rope. 

Figure 2. The evolution of the plasma density and magnetic field in a thermodynamic 
MHD modcl(see text). The plasma density is depicted in color, and projections of the 
magnetic field lines are overlaid on the density. (Corresponding plots of the temperature 
can lie seen in Figure 3.) A closeup of the lower part of the simulation domain near the 
equator is shown, (a) The helmet streamer configuration at the end of the shearing phase, 
(b) The streamer after the magnetic flux has been reduced (reduction of 3.75^ of the 
initial interior bipole). A low-lying filament structure has formed and is just discernible 
near the lower boundary (see Figure 4 for a zoomed view), (e) The filament is now at a 
height of 140. (X)0 km and is moving upward slowly. The enhanced density can he seen to 
he near the bottom of the detached flux surfaces. Flux reduction (at 11.25%) is beyond 
the critical threshold for eruption (sec text), (d) the eruptive phase has started, and 
dense material is carried upward into the corona, shown in (e) and (f). 

Figure 3. The same as Figure 2. except that the evolution of the plasma temperature 
and magnetic field arc' depicted. The high density material shown in Figure 2 is also very 
cold. 


Figuie 4. A zoomed view of the of the prominence-like structure when flux reduction 
is halted at a level of 3.75% (/ - ()08r.,. Figures 2b and 31>) and the calculation is 
continued (see text), allowing a stable flux rope to form. The plasma density (a) and 
temperature (b) are shown together with projections of the magnetic field lines as in 
Figure's 2-3. I racings of the magnetic field lines, colored bv temperature (<■) and density 
(d) are shown. The helical field line's support c old (T ~ 2 x iff 1 A') and dense ( n - 10 KI ) 
material against gravity. 

Figure 5. I he prominence height as function of time for the eruptive phase* of the 
calculation (frames b-f of Figures 2-3). 
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